IMU/UWB Fusion Method Using a Complementary Filter and a Kalman Filter for Hybrid Upper Limb Motion Estimation

Motion capture systems have enormously benefited the research into human–computer interaction in the aerospace field. Given the high cost and susceptibility to lighting conditions of optical motion capture systems, as well as considering the drift in IMU sensors, this paper utilizes a fusion approach with low-cost wearable sensors for hybrid upper limb motion tracking. We propose a novel algorithm that combines the fourth-order Runge–Kutta (RK4) Madgwick complementary orientation filter and the Kalman filter for motion estimation through the data fusion of an inertial measurement unit (IMU) and an ultrawideband (UWB). The Madgwick RK4 orientation filter is used to compensate gyroscope drift through the optimal fusion of a magnetic, angular rate, and gravity (MARG) system, without requiring knowledge of noise distribution for implementation. Then, considering the error distribution provided by the UWB system, we employ a Kalman filter to estimate and fuse the UWB measurements to further reduce the drift error. Adopting the cube distribution of four anchors, the drift-free position obtained by the UWB localization Kalman filter is used to fuse the position calculated by IMU. The proposed algorithm has been tested by various movements and has demonstrated an average decrease in the RMSE of 1.2 cm from the IMU method to IMU/UWB fusion method. The experimental results represent the high feasibility and stability of our proposed algorithm for accurately tracking the movements of human upper limbs.


Introduction
Motion capture (MoCap) systems provide technical support for the operation of space robots, which has enormously benefited the research into human-computer Interaction in the aerospace field. By capturing human upper limb movement and mapping it to space robot motion, space operators can remotely control space robot arms and allow them to perform complex and elaborate actions, such as grasping and handling. Thus, we need to accurately estimate the movement of human upper limbs [1][2][3]. Optical Mocap systems employ multiple cameras to capture the positions of reflective markers and have been widely used in human motion tracking. Although these systems share precise real-time tracking and strong anti-interference abilities, they are expensive and suffer from light and occlusion problems [4,5]. In order to flexibly track human movements, it is necessary to employ ubiquitous wearable sensors.
With the intensive study and rapid advancement of micro-electro-mechanical systems, wearable inertial measurement units (IMUs) have been extensively adopted for 3-D human motion tracking due to their low cost, portability, and high independence [6]. In many of the proposed IMU Mocap methods, accelerometers and gyroscopes are fused to estimate the body segment orientation [7,8], including the knee joint angle [9,10], arm adduction angle [11], and shoulder and elbow joint angles [12,13]. However, it has been illustrated that there is obvious drift error of the attitude angle calculated by the gyroscope, which contributes to the poor estimation accuracy and stability. Considering that the accelerometers can only compensate for the tilt angles relative to the direction of gravity, it is necessary to introduce magnetometers to provide compensation for the yaw angles, which is known as nine-axis IMU or the magnetic, angular rate, and gravity (MARG) system [14]. Data fusion of three inertial sensors can be achieved by adopting either a complementary filter [14][15][16][17] or a classic Kalman filter [18,19]. The Kalman filter has been widely used as an orientation filter for the majority of the proposed methods. Wu et al. [20][21][22][23] employed a particle filter and an unscented Kalman filter for information fusion and developed a micro-sensor Mocap system to achieve real-time tracking. Roetenberg et al. [24] designed a Kalman filter and improved the orientation estimation through the compensation of magnetic disturbances. Considering that the Kalman filter requires knowledge of noise distribution and is prone to computational load and parameter tuning issues, many researchers tend to use the complementary filter. Fourati et al. [25] proposed a complementary observer based on the quaternion method for motion tracking. However, only adopting inertial sensors fusion for motion estimation will still cause drift, so it is expected that a stable and high-precision technology will be introduced to fuse with IMU.
Considering the flexibility and high sampling rate of IMU, ultrawideband (UWB) technology is attractive for data fusion due to its low cost, portability, and low energy consumption [26,27], but it requires a clear line-of-sight (LOS) channel [28]. Compared with traditional UWB positioning methods consisting of maximum likelihood estimation (MLE), linearized least square estimation (LLSE), and weighted centroid estimation (WCE), the extended Kalman filter (EKF) performs with less computational time and less complexity [29]. Applying UWB and IMU fusion for human motion tracking can not only compensate for the low sampling rate of the UWB system, but reduce the drift of IMU as well. Depending on whether the fusion is based on raw time or location, data fusion methods are divided into tightly coupled [30] and loosely coupled [31,32] categories. Kok et al. [33] provided tightly coupled IMU and UWB fusion using an optimization method that showed good performance in terms of handling outliers, but suffered from clock skew challenges. Therefore, in order to facilitate the calculation, the loosely coupled method is often employed. Zihajehzadeh et al. [32,34,35] proposed a Kalman-filter-based IMU and UWB fusion method without a magnetometer and accurately captured lower body motion under magnetic disturbances. However, when extending the rotation matrices to the upper limb, their algorithm may suffer from the singularity problem. Zhang et al. [36] adopted a Mahony filter and quaternion for foot attitude estimation via IMU and UWB fusion, but the use of acceleration double integration to obtain the position would lead to huge cumulative errors.
In this paper, the human upper limb is taken as our research subject, and a novel IMU/UWB data fusion method is proposed for 3-D motion estimation by applying the Runge-Kutta Madgwick filter, the UWB localization Kalman filter, and the IMU/UWB Kalman filter. On the basis of the established kinematics model of the upper limbs, the quaternion method was employed to calculate the attitude angle to avoid the gimbal lock problem. Our proposed algorithm comprises the following two novel aspects: (1) we combined the Madgwick RK4 complementary orientation filter and Kalman filter for motion tracking, the Madgwick RK4 filter was employed to reduce gyroscope drift without leading to noise distribution, and the Kalman filter was implemented for UWB localization and fusion with known error distribution from the UWB system; and (2) the drift-free position obtained by the UWB localization system was used to fuse the position calculated by IMU for upper limb motion estimation, which enormously reduced the drift caused by the double integration of acceleration. The good experimental results represent the high feasibility and stability of our proposed algorithm.
The rest of the paper is structured as follows. In Section 2, the theoretical fusion method for 3-D upper limb tracking is described. The experimental setup and protocol are demonstrated in Section 3. Then, the experimental results are shown and analyzed in Section 4. Finally, conclusions are provided in Section 5.

Theoretical Method
In this section, the proposed information fusion algorithm for human upper limb motion estimation is described. As shown in Figure 1, we first simplified the human upper limb as a kinematics model with three joints (shoulder, elbow, and wrist) and two segments (upper arm and forearm) [21]. The shoulder joint was set as a fixed point with three degrees of freedom (DoFs) and two DoFs for the elbow and wrist [37]. Taking the right arm as an example, two IMU sensors were arranged on the lateral side above the wrist and the elbow, respectively, as well as two UWB tags. The navigation coordinate frame (n) was set as the reference system, which was consistent with the coordinate system of UWB and MoCap. The body frame (b) was attached to each body segment where the sensors were located. It aligned with the n-frame initially. Neglecting the installation errors, the sensor frame (s) was attached to each IMU and aligned with the b-frame to track the movements of human upper limbs. caused by the double integration of acceleration. The good experimental results represent the high feasibility and stability of our proposed algorithm. The rest of the paper is structured as follows. In Section 2, the theoretical fusion method for 3-D upper limb tracking is described. The experimental setup and protocol are demonstrated in Section 3. Then, the experimental results are shown and analyzed in Section 4. Finally, conclusions are provided in Section 5.

Theoretical Method
In this section, the proposed information fusion algorithm for human upper limb motion estimation is described. As shown in Figure 1, we first simplified the human upper limb as a kinematics model with three joints (shoulder, elbow, and wrist) and two segments (upper arm and forearm) [21]. The shoulder joint was set as a fixed point with three degrees of freedom (DoFs) and two DoFs for the elbow and wrist [37]. Taking the right arm as an example, two IMU sensors were arranged on the lateral side above the wrist and the elbow, respectively, as well as two UWB tags. The navigation coordinate frame (n) was set as the reference system, which was consistent with the coordinate system of UWB and MoCap. The body frame (b) was attached to each body segment where the sensors were located. It aligned with the n-frame initially. Neglecting the installation errors, the sensor frame (s) was attached to each IMU and aligned with the bframe to track the movements of human upper limbs. The overall framework of the proposed algorithm consists of three parts, illustrated in Figure 2, including (1) Quaternion RK4 Based Madgwick Orientation Complementary Filter, (2) UWB localization Kalman filter, and (3) IMU/UWB Fusion Kalman filter. For the purpose of estimating the 3-D spatial trajectory of the movements of human upper limbs, we first employed the quaternion method to calculate the attitude angle using gyroscope measurements. However, the integration of gyroscope measurement errors contributed to an accumulating error in the quaternion algorithm [14]. Therefore, we considered adopting the Madgwick orientation filter to improve the motion estimation accuracy through the optimal fusion of the accelerometer, gyroscope, and magnetometer of the MARG system. Then, in order to further reduce the drift error, we proposed that the UWB localization system be combined with the IMU sensors. An extended Kalman filter was utilized to fuse IMU and UWB in order to perform stable and high-precision motion estimation of the human upper limbs. For the purpose of estimating the 3-D spatial trajectory of the movements of human upper limbs, we first employed the quaternion method to calculate the attitude angle using gyroscope measurements. However, the integration of gyroscope measurement errors contributed to an accumulating error in the quaternion algorithm [14]. Therefore, we considered adopting the Madgwick orientation filter to improve the motion estimation accuracy through the optimal fusion of the accelerometer, gyroscope, and magnetometer of the MARG system. Then, in order to further reduce the drift error, we proposed that the UWB localization system be combined with the IMU sensors. An extended Kalman filter was utilized to fuse IMU and UWB in order to perform stable and high-precision motion estimation of the human upper limbs.   Rotation transformation can be described as a vector rotating around a specified rotation axis with a certain angle in a coordinate system. The unit quaternion is defined as q = q 1 + iq 2 + jq 3 + kq 4 (1) For the same vector n r and b r, defined in the n-frame and b-frame, respectively, n r and b r are their extended forms, containing a 0 inserted as the real part [14]. n b q describes the rotation of the n-frame relative to the b-frame, so n r can be expressed as where q * is the conjugate of q; the symbol ⊗ denotes the quaternion product [21], as defined in (3); M n b q is the left multiplication matrix of the quaternion n b q; and M n b q * is the right multiplication matrix of the conjugate quaternion n b q * . Furthermore, the rotation can be represented by a rotation matrix.
where n b C represents the rotation of the n-frame relative to the b-frame.
The quaternion update algorithm uses the angular velocity increment in the sample period measured by the IMU sensors to calculate the quaternion at each time in order to update the human motion data. The quaternion derivative is given by where b ω = 0 b ω x b ω y b ω z T is the angular velocity measured by the gyroscope.
Given the initial value and the rotational angular velocity, the orientation at each time q ω,t . q ω,t ∆t = n bq est,t−1 + where ∆t is the sampling period, n bq est,t−1 is the estimation of the quaternion at the previous time, and b ω t is the angular velocity at time t.
Although increasing the order of the above algorithm can improve the computational accuracy, the complexity is also increasing. The Runge-Kutta method is a high-precision single-step algorithm; the most classic and widely used is the fourth-order Runge-Kutta algorithm. It can perform iterative operations on definite solution problems with known initial values and equations without solving differential equations, which enormously reduces the computational complexity. In this paper, the calculation of human upper limb motion can be regarded as the initial value problem of a differential equation. The calculation formulas are as follows.
where b ω t−1/2 is the angular velocity at the intermediate time between t − 1 and t. The quaternion fourth-order Runge-Kutta algorithm is designed to interpolate in the integration interval, and the slope is iteratively optimized at each step of the calculation to obtain an updated value.

Madgwick RK4 Orientation Complementary Filter for MARG
By formulating an objective function, this filter fuses the data of the tri-axis accelerometer, gyroscope, and magnetometer of the MARG system, then iteratively optimizes it to calculate the orientation n bq . Let the predefined reference vector in the n-frame be nd and the sensor measurement in the b-frame be bŝ . The objective function is defined as By recording the initial quaternion n bq 0 and step size µ, the estimation of n bq n+1 is obtained for n iterations, adopting the general gradient descent algorithm illustrated in (10). The symbol ∇ represents the gradient of the objective function, which can be computed by the objective function and its Jacobian, as shown in (11).
The specific algorithm expression of the accelerometer is described in (12). Normalized gravity nĝ and normalized accelerometer measurements bâ substituted nd and bŝ , respectively.
The accelerometer compensated the tilt angle of motion estimation obtained by the gyroscope measurement, but the heading angle was not compensated, as this requires the introduction of a magnetometer. The magnetometer can be used to measure the strength and direction of the magnetic field and to determine the orientation of a device.
As for the magnetometer, the geomagnetic field nb , as the substitution of nd , can be decomposed into a horizontal axis and a vertical axis theoretically. Additionally, the normalized magnetometer measurement bm is the substitution of bŝ . The specific form can be written as However, the magnetometer may be disturbed by the bias of hard iron and soft iron, causing errors in the measurement direction of the earth's magnetic field, so magnetic distortion compensation needs to be employed.
where nĥ t is the normalized magnetometer measurement in the n-frame and nb t is the compensated geomagnetic field at time t.
Combining the specific algorithms of the accelerometer and magnetometer, the formulas were developed as follows.
By substituting (16) into (11), the gradient of the combined objective function can be written as And the estimated quaternion n b q ∇,t calculated at time t can be given by Employing the Madgwick orientation filter, an estimated rotation n b q est,t was obtained by fusing n b q ω,t and n b q ∇,t .
where γ t represents weight, ensuring that the weighted divergence of n b q ω,t equals the weighted convergence of n b q ∇,t .
According to [14], it is known that when noise augmentation of µ t is assumed to be extremely high, (19) can be rewritten as q ω,t − β ∇f ∇f (20) where β is the divergence rate of n b q ω,t . By substituting (8) into (20), it can be written in the fourth-order Runge-Kuta form as n b q est,t = n bq est,t−1 + ∆t The process of the quaternion RK4-based Madgwick orientation filter algorithm is summarized as follows in Algorithm 1.

UWB Localization Kalman Filter
The UWB positioning system contained a tag-anchor wireless communication channel, as shown in Figure 3. The calculation of the distances between each tag and anchor exploiting double-sided two-way ranging (DS-TWR) method and position estimation employed distances [28,29]. In addition to classical UWB positioning algorithms such as MLE, LLSE, and WCE, we used EKF for position estimation due to its lower complexity and shorter computational time [29].  It was defined that the UWB localization system was aligned with the n-frame. The system model adopted for UWB EKF algorithm can be expressed as It was defined that the UWB localization system was aligned with the n-frame. The system model adopted for UWB EKF algorithm can be expressed as where the subscript u represents the UWB localization Kalman filter, and the state vector the distances between each anchor and tag. The expansion form of h u,t x u,t is expressed in (23). The process noise w u,t−1 and measurement noise v u,t are zero mean additional Gaussian white noise with covariance matrices of Q u,t−1 and R u,t respectively. These matrices were obtained from the UWB system, and distinguished the ranging accuracy of different anchors.
where x i , y i , z i , i = 1, 2, 3, 4 represent the position of the i-th anchor. Thus, the Jacobian matrix H u,t could be calculated as shown in (24). And for the sake of brevity, D i,t , i = 1, 2, 3, 4 is the corresponding row of h u,t x u,t .
Assuming that P u,0 is the initial state estimation covariance and is known as the prior, as well as the initial, state vector x u,0 , the posterior estimations P u,t and x u,t were recursively obtained by employing the prediction and update functions. It was defined that P u,t|t−1 and x u,t|t−1 were the prediction forms illustrated in (25) and could be calculated using the posterior estimations P u,t−1 and x u,t−1 at time t − 1.
Then, the prediction and the measurement vector y u,t were used to update the prior estimations, and the posterior estimations at time t could be expressed as where x u,t is the 3-D position of the UWB tag, and will be used as the measurement vector in the IMU/UWB Kalman filter.

IMU/UWB Fusion Kalman Filter
The innovative aspect of this IMU/UWB Kalman filter is that it uses the drift-free position calculated by the UWB system to compensate for the orientation and position estimated by the IMU system. The information obtained by the IMU system was employed as the state vector, and the position calculated by the UWB system was applied as the measurement vector. The system model can be expressed as where the subscript iu represents the IMU/UWB Kalman filter. x iu,t is the state vector, u iu,t−1 is the input vector, y iu,t is the measurement vector, and w iu,t−1 and v iu,t are the process noise and measurement noise, respectively. The state vector x iu,t = P I MU,t n b q est,t T is a 7 × 1 vector where P I MU,t represents the 3-D position of human upper limb and n b q est,t is the orientation obtained from the Madgwick complementary filter. The state model consists of two parts, and the part n b q est,t can be expressed as where b ω b,t is the bias of the gyroscope and ∇f ∇f is the input vector. Supposing that P I MU,0 is the initial position of the IMU sensors, then P I MU,t can be recursively calculated from a geometric perspective using the rotation quaternion n b q est,t or the rotation matrix n b C t in (5), updated by n b q est,t . The position part of the state model is given by (29) and (30).
where P u I MU,t and n b C u,t represent the position and rotation matrix of the upper arm. The position of the forearm relies on the upper arm, and is given by where P f I MU,t and n b C f ,t represent the position and rotation matrix of the forearm. When using the Kalman filter for the forearm, the input vector should contain the upper arm position.
The Jacobian matrices F iu,t−1 and L iu,t−1 were obtained.
The measurement vector y iu,t = [P UWB,t ] T is a 3 × 1 vector obtained from the UWB localization Kalman filter, representing the 3-D position of UWB tags in the n-frame. The measurement equation can be expressed as where H iu,t is the Jacobian matrix, expressed as The covariance matrices Q iu,t−1 and R iu,t of the process noise w iu,t−1 and the measurement noise v iu,t were calculated as It was assumed that P iu,0 and x iu,0 were the initial state estimation covariance and state vector, respectively, and the prediction forms P iu,t|t−1 and x iu,t|t−1 could be recursively calculated using the following equations.
Then, the prior estimations were updated, and the posterior estimations P iu,t and x iu,t were written as where x iu,t is the fused 3-D position in the IMU/UWB Kalman filter, which was used for 3-D motion reconstruction of human upper limbs. A detailed flowchart of the UWB localization and the IMU/UWB Kalman filters is shown in Figure 4. The parameters were set as follows.
The initial state estimation covariance P u,0 was set to 0.1 × I 3×3 . All elements of P iu,0 were set to 0.1. x u,0 and x iu,0 were obtained from the UWB system and the optical MoCap system for each trial of movement, respectively.
It was assumed that Then, the prior estimations were updated, and the posterior estimations , 0 where , iu t x is the fused 3-D position in the IMU/UWB Kalman filter, which was used for Figure 4. The parameters were set as follows. The initial state estimation covariance ,0 u P was set to 3 3 0.1   I . All elements of ,0 iu P were set to 0.1.

Experimental Setup
To estimate the movement of human upper limbs, wearable 9-axis MPU9250 inertial sensor units (each sensor unit included a tri-axis accelerometer, a tri-axis gyroscope, and a tri-axis magnetometer) were adopted in our experiments. To avoid the extrusive effect of muscles of the human upper limbs, these two IMU sensors were placed on the lateral side, above the wrist and the elbow, respectively [11,21]. Each IMU sensor was connected to the computer via a serial port module for data transmission. The sample rate of the IMU sensors was set to 100 Hz, and the installation direction is shown in Figure 1.
Additionally, a DW1000 UWB real-time localization system manufactured by Haoru Technology, Dalian, China, was employed to track the 3-D position with an update rate of about 10 Hz. The UWB communication technology was implemented between the anchors and tags through the DS-TWR method based on time-domain transmission of radio signals, and the anchors were also connected to the computer via a serial port module. It provided an accuracy of 10 cm for the X/Y axis and 30 cm for the Z axis. Two UWB tags were fixed on the bracelets with the IMU sensors to avoid relative movement, as shown in Figure 5. To ensure the positioning stability of UWB, four anchors were arranged in our laboratory in a cube distribution.

Experimental Setup
To estimate the movement of human upper limbs, wearable 9-axis MPU9250 inertial sensor units (each sensor unit included a tri-axis accelerometer, a tri-axis gyroscope, and a tri-axis magnetometer) were adopted in our experiments. To avoid the extrusive effect of muscles of the human upper limbs, these two IMU sensors were placed on the lateral side, above the wrist and the elbow, respectively [11,21]. Each IMU sensor was connected to the computer via a serial port module for data transmission. The sample rate of the IMU sensors was set to 100 Hz, and the installation direction is shown in Figure 1.
Additionally, a DW1000 UWB real-time localization system manufactured by Haoru Technology, Dalian, China, was employed to track the 3-D position with an update rate of about 10 Hz. The UWB communication technology was implemented between the anchors and tags through the DS-TWR method based on time-domain transmission of radio signals, and the anchors were also connected to the computer via a serial port module. It provided an accuracy of 10 cm for the X/Y axis and 30 cm for the Z axis. Two UWB tags were fixed on the bracelets with the IMU sensors to avoid relative movement, as shown in Figure 5. To ensure the positioning stability of UWB, four anchors were arranged in our laboratory in a cube distribution. The Nokov optical MoCap system captures human motion through sixteen cameras evenly arranged in the laboratory, as illustrated in Figure 5, and is used as a reference system for algorithm verification and further comparison [34]. For each segment of the human upper limb, at least four reflective markers were attached to the skin surface and sensors, constituting an envelope rather than a cluster [9,11,38,39]. The optical MoCap system achieved sub-millimeter positioning accuracy and higher reliability.

Experimental Protocol
Our experiments aimed to capture several common movements of the human upper limbs, including simple whole-arm movements and combined upper-arm and forearm movements. These common movements comprised motion 1: shoulder flexion/extension and abduction/adduction (the whole arm moving along a circular trajectory); motion 2: shoulder abduction/adduction, internal/external rotation, and elbow flexion/extension; and motion 3: shoulder flexion/extension, internal/external rotation, and elbow flexion/extension, as shown in Figure 6. In each trial, one type of movement was performed periodically, with a duration of about one minute. The Nokov optical MoCap system captures human motion through sixteen cameras evenly arranged in the laboratory, as illustrated in Figure 5, and is used as a reference system for algorithm verification and further comparison [34]. For each segment of the human upper limb, at least four reflective markers were attached to the skin surface and sensors, constituting an envelope rather than a cluster [9,11,38,39]. The optical MoCap system achieved sub-millimeter positioning accuracy and higher reliability.

Experimental Protocol
Our experiments aimed to capture several common movements of the human upper limbs, including simple whole-arm movements and combined upper-arm and forearm movements. These common movements comprised motion 1: shoulder flexion/extension and abduction/adduction (the whole arm moving along a circular trajectory); motion 2: shoulder abduction/adduction, internal/external rotation, and elbow flexion/extension; and motion 3: shoulder flexion/extension, internal/external rotation, and elbow flexion/extension, as shown in Figure 6. In each trial, one type of movement was performed periodically, with a duration of about one minute. As for the UWB system, four anchors arranged in conventional cube distribution were employed, as shown in Figure 7. The above movements were performed by one subject under each distribution throughout the experiment. The subject was a 24-year-old healthy female 165 cm in height and 40 kg in weight. It is important to mention that the optical MoCap system tracked the movement of the human upper limbs for each trial, which could be used to initialize the IMU sensors and obtain the positions of UWB anchors. Before each trial, the IMU sensors were calibrated, and a few seconds of stillness then enhanced the stability of the proposed algorithm.

Performance of Madgwick RK4 Orientation Filter
The performance of the Madgwick orientation filter in terms of attitude angle and 3-D position estimation is shown in this section. Figure 8 shows the Euler angles of the wrist for motion 1, the whole arm moving along a circular trajectory, which is beneficial for us in verifying the tri-axis rotation (roll, pitch, and yaw) of the upper limbs. On the basis of the attitude angle obtained by the quaternion fourth-order Runge-Kutta algorithm and the original second-order update algorithm illustrated in (7), the Madgwick orientation filter was combined, as shown in (21) and (20), respectively, and the tri-axis positions of sensor unit 1 during the movement were calculated as shown in Figure 9. The solid lines represent the positions calculated using various methods from the IMU system, whereas the black dotted lines are the reference positions obtained by the optical Mocap system. As can be seen, the proposed methods were able to track the movements of the upper limbs in spite of the errors at the inflection points. As for the UWB system, four anchors arranged in conventional cube distribution were employed, as shown in Figure 7. The above movements were performed by one subject under each distribution throughout the experiment. The subject was a 24-year-old healthy female 165 cm in height and 40 kg in weight. It is important to mention that the optical MoCap system tracked the movement of the human upper limbs for each trial, which could be used to initialize the IMU sensors and obtain the positions of UWB anchors. Before each trial, the IMU sensors were calibrated, and a few seconds of stillness then enhanced the stability of the proposed algorithm. As for the UWB system, four anchors arranged in conventional cube distribution were employed, as shown in Figure 7. The above movements were performed by one subject under each distribution throughout the experiment. The subject was a 24-year-old healthy female 165 cm in height and 40 kg in weight. It is important to mention that the optical MoCap system tracked the movement of the human upper limbs for each trial, which could be used to initialize the IMU sensors and obtain the positions of UWB anchors. Before each trial, the IMU sensors were calibrated, and a few seconds of stillness then enhanced the stability of the proposed algorithm.

Performance of Madgwick RK4 Orientation Filter
The performance of the Madgwick orientation filter in terms of attitude angle and 3-D position estimation is shown in this section. Figure 8 shows the Euler angles of the wrist for motion 1, the whole arm moving along a circular trajectory, which is beneficial for us in verifying the tri-axis rotation (roll, pitch, and yaw) of the upper limbs. On the basis of the attitude angle obtained by the quaternion fourth-order Runge-Kutta algorithm and the original second-order update algorithm illustrated in (7), the Madgwick orientation filter was combined, as shown in (21) and (20), respectively, and the tri-axis positions of sensor unit 1 during the movement were calculated as shown in Figure 9. The solid lines represent the positions calculated using various methods from the IMU system, whereas the black dotted lines are the reference positions obtained by the optical Mocap system. As can be seen, the proposed methods were able to track the movements of the upper limbs in spite of the errors at the inflection points.

Performance of Madgwick RK4 Orientation Filter
The performance of the Madgwick orientation filter in terms of attitude angle and 3-D position estimation is shown in this section. Figure 8 shows the Euler angles of the wrist for motion 1, the whole arm moving along a circular trajectory, which is beneficial for us in verifying the tri-axis rotation (roll, pitch, and yaw) of the upper limbs. On the basis of the attitude angle obtained by the quaternion fourth-order Runge-Kutta algorithm and the original second-order update algorithm illustrated in (7), the Madgwick orientation filter was combined, as shown in (21) and (20), respectively, and the tri-axis positions of sensor unit 1 during the movement were calculated as shown in Figure 9. The solid lines represent the positions calculated using various methods from the IMU system, whereas the black dotted lines are the reference positions obtained by the optical Mocap system. As can be seen, the proposed methods were able to track the movements of the upper limbs in spite of the errors at the inflection points.  The position accuracies of the three kinds of movements are summarized in Table 1. It shows the 3-D position accuracy of forearm motion tracking by the fourth-order Runge-Kutta Madgwick method and the original second-order update algorithm compared with the optical Mocap system. According to Table 1, the root mean square error (RMSE) values of tri-axis motion tracking of Madgwick RK4 method varied from 6.17 cm to 12.76 cm, which shows the feasibility of precise motion tracking of the human upper limbs. In order to further reduce the error, it is necessary to utilize the UWB measurements.   The position accuracies of the three kinds of movements are summarized in Table 1. It shows the 3-D position accuracy of forearm motion tracking by the fourth-order Runge-Kutta Madgwick method and the original second-order update algorithm compared with the optical Mocap system. According to Table 1, the root mean square error (RMSE) values of tri-axis motion tracking of Madgwick RK4 method varied from 6.17 cm to 12.76 cm, which shows the feasibility of precise motion tracking of the human upper limbs. In order to further reduce the error, it is necessary to utilize the UWB measurements.  The position accuracies of the three kinds of movements are summarized in Table 1. It shows the 3-D position accuracy of forearm motion tracking by the fourth-order Runge-Kutta Madgwick method and the original second-order update algorithm compared with the optical Mocap system. According to Table 1, the root mean square error (RMSE) values of tri-axis motion tracking of Madgwick RK4 method varied from 6.17 cm to 12.76 cm, which shows the feasibility of precise motion tracking of the human upper limbs. In order to further reduce the error, it is necessary to utilize the UWB measurements.

Performance of UWB Localization Kalman Filter
As for the UWB localization Kalman filter, we utilized multiple sets of movements of the human upper limb to verify the accuracy of the UWB positioning system. It is important to mention that the installation direction of the IMU was always the same as that in Figure 1 during the experiment. And the subject stood at the edge of the UWB area to ensure that there was no occlusion between the tags and anchors during the movement of the upper limbs. Before the experiment, we first kept the tag still and verified the static positioning accuracy of UWB. Table 2 shows the static positioning variance of the four-anchor UWB system. The mean values of the tri-axis variances fell within 0.6 cm 2 , so it shared the good positioning robustness of the UWB system.  Figure 10 shows a comparison of the position of upper limb motion tracking for three kinds of movements by the IMU, UWB, and optical Mocap systems. Motion 1, represented by Figure 10a, consisted of movement along a circular trajectory, which is consistent with of the motion in Figure 9. And motion 2 was a combined motion, represented by Figure 10b, comprising shoulder abduction/adduction, internal/external rotation, and elbow flexion/extension. Motion 3 was another combined motion, represented by Figure 10c and similar to motion 2, consisting of shoulder flexion/extension, internal/external rotation, and elbow flexion/extension. As can be seen in these figures, it is obvious that compared with the position calculated by the IMU sensors, the position calculated by the UWB system was closer to the reference value of the optical Mocap system at certain steps, especially the inflection points. Therefore, as for the measurement, the UWB system can be fused with IMU for upper limb motion estimation to reduce the error caused by the drift of the IMU sensors.

Performance of UWB Localization Kalman Filter
As for the UWB localization Kalman filter, we utilized multiple sets of movements of the human upper limb to verify the accuracy of the UWB positioning system. It is important to mention that the installation direction of the IMU was always the same as that in Figure 1 during the experiment. And the subject stood at the edge of the UWB area to ensure that there was no occlusion between the tags and anchors during the movement of the upper limbs. Before the experiment, we first kept the tag still and verified the static positioning accuracy of UWB. Table 2 shows the static positioning variance of the fouranchor UWB system. The mean values of the tri-axis variances fell within 0.6 cm 2 , so it shared the good positioning robustness of the UWB system. Figure 10 shows a comparison of the position of upper limb motion tracking for three kinds of movements by the IMU, UWB, and optical Mocap systems. Motion 1, represented by Figure 10a, consisted of movement along a circular trajectory, which is consistent with of the motion in Figure 9. And motion 2 was a combined motion, represented by Figure  10b, comprising shoulder abduction/adduction, internal/external rotation, and elbow flexion/extension. Motion 3 was another combined motion, represented by Figure 10c and similar to motion 2, consisting of shoulder flexion/extension, internal/external rotation, and elbow flexion/extension. As can be seen in these figures, it is obvious that compared with the position calculated by the IMU sensors, the position calculated by the UWB system was closer to the reference value of the optical Mocap system at certain steps, especially the inflection points. Therefore, as for the measurement, the UWB system can be fused with IMU for upper limb motion estimation to reduce the error caused by the drift of the IMU sensors.

Performance of IMU/UWB Kalman Filter
The performance of the IMU/UWB Kalman filter in 3-D upper limb motion tracking is shown in this section. As mentioned previously, by adopting the cube distribution of four anchors, the measurement of UWB was fused with IMU for drift error reduction and better positioning accuracy in motion estimation. The position's RMSE values are illustrated in Table 3, where motion 1, motion 2, and motion 3 are represented in Figure 10a-c, respectively. This shows a comparison of the 3-D position accuracy of various upper limb motion tracking methods, i.e., the IMU method, UWB method, and IMU/UWB fusion method with the optical Mocap system. Although the position accuracy of the UWB system may have sometimes fallen lower than that of the IMU, the adopted Kalman filter adaptively adjusted the weights between them based on the covariances. It can be observed that the tri-axis position RMSE calculated by the IMU/UWB Kalman filter was consistently less than that without UWB fusion. In our experiment, the proposed method was been extensively tested by various movements. For simple movements like motion 1, consisting of a circular trajectory, the tri-axis position RMSE values were all less than 10 cm. Compared with the IMU method, the relative errors calculated by the IMU/UWB fusion method were reduced by 40%, 3.6%, and 25.5% in the X-axis, Y-axis, and Z-axis, respectively. Complex combined movements such as motion 2 and 3 always consisted of multiple simple movements, such as shoulder flexion/extension, abduction/adduction, internal/external rotation, and elbow flexion/extension. Despite the fact that the tri-axis errors with combined motions are usually more significant than those with simple movements, the IMU/UWB fusion method still demonstrated better position accuracy than the IMU method, with a maximum RMSE of 12.2 cm. And our proposed methodology achieved an average decrease in the RMSE of 1.2 cm from the IMU method to the IMU/UWB fusion method. In comparison, by transferring the angle RMSE into position RMSE, the position RMSE values illustrated in Table 3 fell within the accuracy range of forearm motion without the constraints represented in [21].  Figure 11 shows a comparison of the 3-D spatial trajectory of upper limb motion reconstructed by the proposed algorithm and the real movement for motion 3. As can be seen in this figure, the proposed method was able to accurately track the movement of the human upper limb, and showed high feasibility and stability.

Power Consumption and Cost
According to the datasheet of the MPU9250 which was adopted in our experiment, the power consumption of the sensor depends on the selected operating mode and sensor configuration. In general, the device consumes between 8.3 μA and 3.9 mA. When the voltage is 3.3 V, the power consumption is as follows: in low-power mode, it ranges from 27.4 μW to 66.7 μW; in normal mode, with both the accelerometer and the gyroscope enabled, it varies from 31.5 μW to 764.7 μW; and in normal mode with all three sensors (accelerometer, gyroscope, and magnetometer) enabled, the power consumption ranges

Power Consumption and Cost
According to the datasheet of the MPU9250 which was adopted in our experiment, the power consumption of the sensor depends on the selected operating mode and sensor configuration. In general, the device consumes between 8.3 µA and 3.9 mA. When the voltage is 3.3 V, the power consumption is as follows: in low-power mode, it ranges from 27.4 µW to 66.7 µW; in normal mode, with both the accelerometer and the gyroscope enabled, it varies from 31.5 µW to 764.7 µW; and in normal mode with all three sensors (accelerometer, gyroscope, and magnetometer) enabled, the power consumption ranges from 37.3 µW to 926.7 µW. The average price for a single MPU9250 chip is between USD 1 and 5.
As for the UWB system, according to the DW1000 datasheet, the typical power consumption values for the device in different modes are as follows: in idle mode, 24.75 mW; in receive mode, 56.84 mW; in transmit mode (at 6.8 Mbps data rate), around 270.76 mW; and in sleep mode, 20 nW. Typically, the cost of a single DW1000 module ranges from around USD 10 to 25.
However, the power consumption of an optical MoCap system can range from a few hundred watts to several kilowatts. This is because the cameras employed in this system require significant processing power and strict light conditions. And the cost of a basic system with a few cameras and basic software is around USD 1000. Therefore, compared with the optical MoCap system, the IMU/UWB system which we adopted exhibits extremely lower power consumption and is more economical and applicable.

Conclusions
In this paper, we proposed a novel method for hybrid upper limb motion tracking and 3-D positioning by fusing the IMU and UWB systems. We first simplified the human upper limb as a kinematics model and employed the quaternion method to calculate the attitude angle of each segment. To compensate for the accumulated error of gyroscope measurement, the fourth-order Runge-Kutta Madgwick orientation filter was adopted to improve the accuracy of 3-D motion tracking through the optimal fusion of the accelerometer, gyroscope, and magnetometer of the MARG system. The RMSE values varied from 6.17 cm to 12.76 cm, which shows feasibility for the precise motion tracking of human upper limbs. The error was mainly caused by drift of the IMU sensors. And it was inevitable that there would be a slight misalignment of coordinate frames at the initial moment.
In order to further reduce the drift error, we combined the UWB localization system with the IMU sensors. The static positioning variance of the four-anchor UWB system was tested, and the mean values of the tri-axis variances were within 0.6 cm 2 , so it shared good positioning robustness. By employing the UWB localization Kalman filter, the accuracy of UWB was verified by multiple sets of movements of the human upper limbs.
Adopting the four-anchor UWB system, we employed various movements to test the IMU/UWB Kalman filter. The experimental results represent that our proposed fusion algorithm achieved an average decrease in the RMSE of 1.2 cm from the IMU method to the IMU/UWB fusion method. With high feasibility and stability, our proposed algorithm was able to accurately track the movements of the human upper limbs.
In future work, we aim to study the effects of various distributions of multiple anchors on UWB localization accuracy, as well as to extend the proposed algorithm to the wholebody motion estimation.